home *** CD-ROM | disk | FTP | other *** search
/ Celestin Apprentice 2 / Apprentice-Release2.iso / Source Code / C / Applications / RTrace 1.0 / source / noise.c < prev    next >
Encoding:
C/C++ Source or Header  |  1992-08-16  |  11.0 KB  |  400 lines  |  [TEXT/KAHL]

  1. /*
  2.  * Copyright (c) 1988, 1992 Antonio Costa, INESC-Norte.
  3.  * All rights reserved.
  4.  *
  5.  * This code received contributions from the following people:
  6.  *
  7.  *  Roman Kuchkuda      - basic ray tracer
  8.  *  Mark VandeWettering - MTV ray tracer
  9.  *  Augusto Sousa       - overall, shading model
  10.  *  Craig Kolb          - noise code
  11.  *  Robert Skinner      - noise functions
  12.  *  David Buck          - noise functions
  13.  *  Reid Judd        - corrections
  14.  *
  15.  * Redistribution and use in source and binary forms are permitted
  16.  * provided that the above copyright notice and this paragraph are
  17.  * duplicated in all such forms and that any documentation,
  18.  * advertising materials, and other materials related to such
  19.  * distribution and use acknowledge that the software was developed
  20.  * by Antonio Costa, at INESC-Norte. The name of the author and
  21.  * INESC-Norte may not be used to endorse or promote products derived
  22.  * from this software without specific prior written permission.
  23.  * THIS SOFTWARE IS PROVIDED ``AS IS'' AND WITHOUT ANY EXPRESS OR
  24.  * IMPLIED WARRANTIES, INCLUDING, WITHOUT LIMITATION, THE IMPLIED
  25.  * WARRANTIES OF MERCHANTIBILITY AND FITNESS FOR A PARTICULAR PURPOSE.
  26.  */
  27. #include "defs.h"
  28. #include "extern.h"
  29.  
  30. /**********************************************************************
  31.  *    RAY TRACING - Noise - Version 7.3.1                             *
  32.  *                                                                    *
  33.  *    MADE BY    : Antonio Costa, INESC-Norte, March 1990             *
  34.  *    MODIFIED BY: Antonio Costa, INESC-Norte, July 1992              *
  35.  **********************************************************************/
  36.  
  37. #define S_CURVE(a) ((a) * (a) * (3.0 - 2.0 * (a)))
  38.  
  39. #define HASH(a, b, c)\
  40. hash_table[hash_table[hash_table[(a) BIT_AND 4095] BIT_XOR\
  41. ((b) BIT_AND 4095)] BIT_XOR\
  42. ((c) BIT_AND 4095)]
  43.  
  44. #define INCRSUM(m, s, x, y, z)\
  45. ((s) * (random_table[m] * 0.5\
  46. + random_table[m + 1] * (x)\
  47. + random_table[m + 2] * (y)\
  48. + random_table[m + 3] * (z)))
  49.  
  50. #define RANDOM_TABLE_SIZE (267)
  51. #define SIN_TABLE_SIZE (1000)
  52.  
  53. static real     random_table[RANDOM_TABLE_SIZE];
  54. static short int *hash_table;
  55. static real     sin_table[SIN_TABLE_SIZE];
  56. static unsigned short int crc_table[256] =
  57. {
  58.   0x0000, 0xc0c1, 0xc181, 0x0140, 0xc301, 0x03c0, 0x0280, 0xc241,
  59.   0xc601, 0x06c0, 0x0780, 0xc741, 0x0500, 0xc5c1, 0xc481, 0x0440,
  60.   0xcc01, 0x0cc0, 0x0d80, 0xcd41, 0x0f00, 0xcfc1, 0xce81, 0x0e40,
  61.   0x0a00, 0xcac1, 0xcb81, 0x0b40, 0xc901, 0x09c0, 0x0880, 0xc841,
  62.   0xd801, 0x18c0, 0x1980, 0xd941, 0x1b00, 0xdbc1, 0xda81, 0x1a40,
  63.   0x1e00, 0xdec1, 0xdf81, 0x1f40, 0xdd01, 0x1dc0, 0x1c80, 0xdc41,
  64.   0x1400, 0xd4c1, 0xd581, 0x1540, 0xd701, 0x17c0, 0x1680, 0xd641,
  65.   0xd201, 0x12c0, 0x1380, 0xd341, 0x1100, 0xd1c1, 0xd081, 0x1040,
  66.   0xf001, 0x30c0, 0x3180, 0xf141, 0x3300, 0xf3c1, 0xf281, 0x3240,
  67.   0x3600, 0xf6c1, 0xf781, 0x3740, 0xf501, 0x35c0, 0x3480, 0xf441,
  68.   0x3c00, 0xfcc1, 0xfd81, 0x3d40, 0xff01, 0x3fc0, 0x3e80, 0xfe41,
  69.   0xfa01, 0x3ac0, 0x3b80, 0xfb41, 0x3900, 0xf9c1, 0xf881, 0x3840,
  70.   0x2800, 0xe8c1, 0xe981, 0x2940, 0xeb01, 0x2bc0, 0x2a80, 0xea41,
  71.   0xee01, 0x2ec0, 0x2f80, 0xef41, 0x2d00, 0xedc1, 0xec81, 0x2c40,
  72.   0xe401, 0x24c0, 0x2580, 0xe541, 0x2700, 0xe7c1, 0xe681, 0x2640,
  73.   0x2200, 0xe2c1, 0xe381, 0x2340, 0xe101, 0x21c0, 0x2080, 0xe041,
  74.   0xa001, 0x60c0, 0x6180, 0xa141, 0x6300, 0xa3c1, 0xa281, 0x6240,
  75.   0x6600, 0xa6c1, 0xa781, 0x6740, 0xa501, 0x65c0, 0x6480, 0xa441,
  76.   0x6c00, 0xacc1, 0xad81, 0x6d40, 0xaf01, 0x6fc0, 0x6e80, 0xae41,
  77.   0xaa01, 0x6ac0, 0x6b80, 0xab41, 0x6900, 0xa9c1, 0xa881, 0x6840,
  78.   0x7800, 0xb8c1, 0xb981, 0x7940, 0xbb01, 0x7bc0, 0x7a80, 0xba41,
  79.   0xbe01, 0x7ec0, 0x7f80, 0xbf41, 0x7d00, 0xbdc1, 0xbc81, 0x7c40,
  80.   0xb401, 0x74c0, 0x7580, 0xb541, 0x7700, 0xb7c1, 0xb681, 0x7640,
  81.   0x7200, 0xb2c1, 0xb381, 0x7340, 0xb101, 0x71c0, 0x7080, 0xb041,
  82.   0x5000, 0x90c1, 0x9181, 0x5140, 0x9301, 0x53c0, 0x5280, 0x9241,
  83.   0x9601, 0x56c0, 0x5780, 0x9741, 0x5500, 0x95c1, 0x9481, 0x5440,
  84.   0x9c01, 0x5cc0, 0x5d80, 0x9d41, 0x5f00, 0x9fc1, 0x9e81, 0x5e40,
  85.   0x5a00, 0x9ac1, 0x9b81, 0x5b40, 0x9901, 0x59c0, 0x5880, 0x9841,
  86.   0x8801, 0x48c0, 0x4980, 0x8941, 0x4b00, 0x8bc1, 0x8a81, 0x4a40,
  87.   0x4e00, 0x8ec1, 0x8f81, 0x4f40, 0x8d01, 0x4dc0, 0x4c80, 0x8c41,
  88.   0x4400, 0x84c1, 0x8581, 0x4540, 0x8701, 0x47c0, 0x4680, 0x8641,
  89.   0x8201, 0x42c0, 0x4380, 0x8341, 0x4100, 0x81c1, 0x8081, 0x4040
  90. };
  91. static void
  92. init_hash_table()
  93. {
  94.   REG int         i, j;
  95.   short int       temp;
  96.  
  97.   ALLOCATE(hash_table, short int, 4096);
  98.   for (i = 0; i < 4096; POSINC(i))
  99.     hash_table[i] = i;
  100.   for (i = 4095; i > 0; POSDEC(i))
  101.   {
  102.     j = (int) (RANDOM * 4096.0);
  103.     temp = hash_table[i];
  104.     hash_table[i] = hash_table[j];
  105.     hash_table[j] = temp;
  106.   }
  107. }
  108. /*
  109.  * ATTENTION: Note that passing a real to crc16 and interpreting it as
  110.  * an array of chars means that machines with different floating-point
  111.  * representation schemes will evaluate black_noise(point) differently.
  112.  */
  113. static unsigned int
  114. crc16(byte, count)
  115.   REG char_ptr    byte;
  116.   REG int         count;
  117. {
  118.   REG unsigned int crc;
  119.  
  120.   crc = 0;
  121.   while (POSDEC(count))
  122.   {
  123.     crc = (crc SHR 8) BIT_XOR
  124.       crc_table[(unsigned char) (crc BIT_XOR *byte)];
  125.     POSINC(byte);
  126.   }
  127.   return crc;
  128. }
  129. static unsigned int
  130. multiply(p)
  131.   xyz_ptr         p;
  132. {
  133.   p->x *= 0.12345;
  134.   p->y *= 0.12345;
  135.   p->z *= 0.12345;
  136.   return crc16((char_ptr) p, sizeof(xyz_struct));
  137. }
  138. void
  139. init_texture()
  140. {
  141.   REG int         i;
  142.   xyz_struct      p;
  143.  
  144.   init_hash_table();
  145.   for (i = 0; i < RANDOM_TABLE_SIZE; POSINC(i))
  146.   {
  147.     p.x = (real) i;
  148.     p.y = (real) i;
  149.     p.z = (real) i;
  150.     random_table[i] = (real) multiply(&p) / 32768.0 - 1.0;
  151.   }
  152.   for (i = 0; i < SIN_TABLE_SIZE; POSINC(i))
  153.     sin_table[i] = SIN((real) i / (real) SIN_TABLE_SIZE * 2.0 * PI);
  154.   p.y = 0.0;
  155.   p.z = 0.0;
  156.   for (i = 0; i < WAVES_MAX; POSINC(i))
  157.   {
  158.     p.x = (real) i;
  159.     black_noise2(&p, &(wave.source[i]));
  160.     NORMALIZE(wave.source[i]);
  161.     wave.frequency[i] = RANDOM;
  162.   }
  163. }
  164. real
  165. black_noise(p)
  166.   xyz_ptr         p;
  167. {
  168.   REG long        ix, iy, iz, jx, jy, jz;
  169.   real            x, y, z, fx, fy, fz;
  170.   real            sx, sy, sz, tx, ty, tz;
  171.   real            sum;
  172.   short int       m;
  173.  
  174.   x = p->x + MAXINT DIV 2;
  175.   y = p->y + MAXINT DIV 2;
  176.   z = p->z + MAXINT DIV 2;
  177.   ix = (long) x;
  178.   iy = (long) y;
  179.   iz = (long) z;
  180.   jx = ix + 1;
  181.   jy = iy + 1;
  182.   jz = iz + 1;
  183.   fx = FRAC(x);
  184.   fy = FRAC(y);
  185.   fz = FRAC(z);
  186.   sx = S_CURVE(fx);
  187.   sy = S_CURVE(fy);
  188.   sz = S_CURVE(fz);
  189.   tx = 1.0 - sx;
  190.   ty = 1.0 - sy;
  191.   tz = 1.0 - sz;
  192.  
  193.   m = HASH(ix, iy, iz) BIT_AND 255;
  194.   sum = INCRSUM(m, tx * ty * tz, fx, fy, fz);
  195.   m = HASH(jx, iy, iz) BIT_AND 255;
  196.   sum += INCRSUM(m, sx * ty * tz, fx - 1.0, fy, fz);
  197.   m = HASH(ix, jy, iz) BIT_AND 255;
  198.   sum += INCRSUM(m, tx * sy * tz, fx, fy - 1.0, fz);
  199.   m = HASH(jx, jy, iz) BIT_AND 255;
  200.   sum += INCRSUM(m, sx * sy * tz, fx - 1.0, fy - 1.0, fz);
  201.   m = HASH(ix, iy, jz) BIT_AND 255;
  202.   sum += INCRSUM(m, tx * ty * sz, fx, fy, fz - 1.0);
  203.   m = HASH(jx, iy, jz) BIT_AND 255;
  204.   sum += INCRSUM(m, sx * ty * sz, fx - 1.0, fy, fz - 1.0);
  205.   m = HASH(ix, jy, jz) BIT_AND 255;
  206.   sum += INCRSUM(m, tx * sy * sz, fx, fy - 1.0, fz - 1.0);
  207.   m = HASH(jx, jy, jz) BIT_AND 255;
  208.   sum += INCRSUM(m, sx * sy * sz, fx - 1.0, fy - 1.0, fz - 1.0);
  209.   sum = MAX(0.0, MIN(1.0, sum + 0.5));
  210.   return sum;
  211. }
  212. void
  213. black_noise2(p, r)
  214.   xyz_ptr         p, r;
  215. {
  216.   REG long        ix, iy, iz, jx, jy, jz;
  217.   real            x, y, z, fx, fy, fz;
  218.   real            sx, sy, sz, tx, ty, tz;
  219.   real            px, py, pz, s;
  220.   short int       m;
  221.  
  222.   x = p->x + MAXINT DIV 2;
  223.   y = p->y + MAXINT DIV 2;
  224.   z = p->z + MAXINT DIV 2;
  225.   ix = (long) x;
  226.   iy = (long) y;
  227.   iz = (long) z;
  228.   jx = ix + 1;
  229.   jy = iy + 1;
  230.   jz = iz + 1;
  231.   fx = FRAC(x);
  232.   fy = FRAC(y);
  233.   fz = FRAC(z);
  234.   sx = S_CURVE(fx);
  235.   sy = S_CURVE(fy);
  236.   sz = S_CURVE(fz);
  237.   tx = 1.0 - sx;
  238.   ty = 1.0 - sy;
  239.   tz = 1.0 - sz;
  240.  
  241.   m = HASH(ix, iy, iz) BIT_AND 255;
  242.   px = fx;
  243.   py = fy;
  244.   pz = fz;
  245.   s = tx * ty * tz;
  246.   r->x = INCRSUM(m, s, px, py, pz);
  247.   r->y = INCRSUM(m + 4, s, px, py, pz);
  248.   r->z = INCRSUM(m + 8, s, px, py, pz);
  249.   m = HASH(jx, iy, iz) BIT_AND 255;
  250.   px = fx - 1.0;
  251.   s = sx * ty * tz;
  252.   r->x += INCRSUM(m, s, px, py, pz);
  253.   r->y += INCRSUM(m + 4, s, px, py, pz);
  254.   r->z += INCRSUM(m + 8, s, px, py, pz);
  255.   m = HASH(jx, jy, iz) BIT_AND 255;
  256.   py = fy - 1.0;
  257.   s = sx * sy * tz;
  258.   r->x += INCRSUM(m, s, px, py, pz);
  259.   r->y += INCRSUM(m + 4, s, px, py, pz);
  260.   r->z += INCRSUM(m + 8, s, px, py, pz);
  261.   m = HASH(ix, jy, iz) BIT_AND 255;
  262.   px = fx;
  263.   s = tx * sy * tz;
  264.   r->x += INCRSUM(m, s, px, py, pz);
  265.   r->y += INCRSUM(m + 4, s, px, py, pz);
  266.   r->z += INCRSUM(m + 8, s, px, py, pz);
  267.   m = HASH(ix, jy, jz) BIT_AND 255;
  268.   pz = fz - 1.0;
  269.   s = tx * sy * sz;
  270.   r->x += INCRSUM(m, s, px, py, pz);
  271.   r->y += INCRSUM(m + 4, s, px, py, pz);
  272.   r->z += INCRSUM(m + 8, s, px, py, pz);
  273.   m = HASH(jx, jy, jz) BIT_AND 255;
  274.   px = fx - 1.0;
  275.   s = sx * sy * sz;
  276.   r->x += INCRSUM(m, s, px, py, pz);
  277.   r->y += INCRSUM(m + 4, s, px, py, pz);
  278.   r->z += INCRSUM(m + 8, s, px, py, pz);
  279.   m = HASH(jx, iy, jz) BIT_AND 255;
  280.   py = fy;
  281.   s = sx * ty * sz;
  282.   r->x += INCRSUM(m, s, px, py, pz);
  283.   r->y += INCRSUM(m + 4, s, px, py, pz);
  284.   r->z += INCRSUM(m + 8, s, px, py, pz);
  285.   m = HASH(ix, iy, jz) BIT_AND 255;
  286.   px = fx;
  287.   s = tx * ty * sz;
  288.   r->x += INCRSUM(m, s, px, py, pz);
  289.   r->y += INCRSUM(m + 4, s, px, py, pz);
  290.   r->z += INCRSUM(m + 8, s, px, py, pz);
  291. }
  292. real
  293. fbm_noise(p, omega, lambda, octaves)
  294.   REG xyz_ptr     p;
  295.   real            omega, lambda;
  296.   int             octaves;
  297. {
  298.   REG int         i;
  299.   REG real        l, o, value;
  300.   xyz_struct      temp;
  301.  
  302.   value = black_noise(p);
  303.   l = lambda;
  304.   o = omega;
  305.   for (i = 2; i <= octaves; POSINC(i))
  306.   {
  307.     temp.x = l * p->x;
  308.     temp.y = l * p->y;
  309.     temp.z = l * p->z;
  310.     value += o * black_noise(&temp);
  311.     if (i < octaves)
  312.     {
  313.       l *= lambda;
  314.       o *= omega;
  315.     }
  316.   }
  317.   return value;
  318. }
  319. void
  320. fbm_noise2(p, omega, lambda, octaves, r)
  321.   REG xyz_ptr     p;
  322.   real            omega, lambda;
  323.   int             octaves;
  324.   xyz_ptr         r;
  325. {
  326.   REG int         i;
  327.   REG real        l, o;
  328.   xyz_struct      temp, value;
  329.  
  330.   black_noise2(p, r);
  331.   l = lambda;
  332.   o = omega;
  333.   for (i = 2; i <= octaves; POSINC(i))
  334.   {
  335.     temp.x = l * p->x;
  336.     temp.y = l * p->y;
  337.     temp.z = l * p->z;
  338.     black_noise2(&temp, &value);
  339.     r->x += o * value.x;
  340.     r->y += o * value.y;
  341.     r->z += o * value.z;
  342.     if (i < octaves)
  343.     {
  344.       l *= lambda;
  345.       o *= omega;
  346.     }
  347.   }
  348. }
  349. real
  350. chaos_noise(p, octaves)
  351.   xyz_ptr         p;
  352.   int             octaves;
  353. {
  354.   return fbm_noise(p, 0.5, 2.0, octaves);
  355. }
  356. real
  357. marble_noise(p)
  358.   xyz_ptr         p;
  359. {
  360.   REG real        t;
  361.  
  362.   t = (SIN(8.0 * chaos_noise(p, 6) + 7.0 * p->z) + 1.0) * 0.5;
  363.   return POWER(t, 0.77);
  364. }
  365. real
  366. turbulence(p)
  367.   xyz_ptr         p;
  368. {
  369.   return fbm_noise(p, 0.5, 2.0, 4);
  370. }
  371. void
  372. turbulence2(p, r)
  373.   xyz_ptr         p, r;
  374. {
  375.   fbm_noise2(p, 0.5, 2.0, 4, r);
  376. }
  377. real
  378. cycloidal(value)
  379.   real            value;
  380. {
  381.   REG real        t;
  382.  
  383.   t = sin_table[(unsigned int) (FRAC(ABS(value)) * (real) SIN_TABLE_SIZE)];
  384.   return t;
  385. }
  386. real
  387. triangle_wave(value)
  388.   real            value;
  389. {
  390.   REG real        offset;
  391.  
  392.   if (value > 0.0)
  393.     offset = FRAC(value);
  394.   else
  395.     offset = 1.0 - FRAC(ABS(value));
  396.   if (offset > 0.5)
  397.     return (2.0 * (1.0 - offset));
  398.   return (2.0 * offset);
  399. }
  400.